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Abstract 

The perturbation theory based on the Riemann-Hilbert problem is developed for the 
modified nonlinear Schrodinger equation which describes the propagation of femtosecond 
optical pulses in nonlinear single-mode optical fibers. A detailed analysis of the adiabatic 
approximation to perturbation-induced evolution of the soliton parameters is given. The 
linear perturbation and the Raman gain are considered as examples. 

PACS. 03.40Kf - Waves and wave propagation: general mathematical aspects. 02.30Jr - 
Partial differential equations. 
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1. Introduction 

The study of the dynamical processes associated with the propagation of high-power optical pulses in 
single-mode nonlinear fibers is based as a rule on the integrable nonlinear Schrodinger equation (NLSE) [1]. 
Various realistic effects accompanying the soliton propagation and destroying the integrability of the NLSE 
are usually treated as perturbations. There are different approaches to describe analytically perturbation- 
induced dynamics of the NLS solitons [2-8], for review see Ref. 9. ft is evident that the "quality" of 
taking into account for the above effects depends crucially on smallness of a parameter responsible for a 
definite perturbation. Just this situation takes place with the Kerr nonlinearity dispersion effect. Being 
sufficiently small in the picosecond pulse duration region, it becomes essential for femtosecond solitons, 
having a parameter of the order f0~ 2 — 10 . Hence, streactly speaking, this effect cannot be treated as a 
perturbation for the femtosecond region of soliton pulse duration. 
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The natural approach to treat analytically the dynamics of femtosecond solitons is to consider the so- 
called perturbed modified nonlinear Schrodingcr equation (MNLSE) [10] 



iq z + 2<7rr + ia (\q\ 2 q) T + \q\ 2 q = r, 



(l) 



where the term with the real parameter a governs the effect of the Kerr nonlinearity dispersion (self- 
steepening) and r accounts for small effects which we will consider as perturbation. Here q(r, z) is the 
normalized slowly varying amplitude of the complex field envelope, z is the normalized propagation distance 
along the fiber, r is the normalized time measured in a frame of reference moving with the pulse at the 
group velocity (the retarded time). It is remarkable that MNLSE (1) with zero r.h.s. is still integrable by 
the inverse scattering transform (1ST) method [11], though the linear spectral problem associated with the 
MNLSE is different from that for the NLSE. 

Our primary goal is to develop a simple formalism to treat analytically the femtosecond soliton dynam- 
ics governed by Eq. (1). Three points should be stressed which differ our approach from the previously 
known ones. First, we account for the Kerr nonlinearity dispersion effect exactly. In other words, we do not 
make any hypothesis about smallness of a in Eq. (1). Moreover, we consider as a background solution not 
the sech-like pulse of the NLS type but precisely the MNLS soliton. Finally, the third point is relevant to 
the formalism, namely, we employ the Riemann-Hilbert (RH) problem, which was proved to be effective for 
treating perturbations to nonlinear evolution equations integrable by means of the Zakharov-Shabat spectral 
problem [12-15]. Recently, we developed the RH problem-based approach [16-18] for solving nonlinear equa- 
tions integrable by the Wadati-Konno-Ichikawa spectral problem [11]. This approach includes the MNLSE 
and some its generalizations. This development serves as a base for taking into account small perturbations. 

The paper is arranged as follows. In Sec. 2 the basic results on the RH problem-based approach to 
the MNLSE is summarized. In Sec. 3 the general one-soliton solution of the MNLSE is derived in a form 
which, we believe, is as simple as possible. Here we also discuss the limiting transition to the NLS soliton. 
In Sec. 4 we obtain the perturbation-induced evolution equations for the RH problem data related to the 
soliton parameters and discuss peculiarities of the perturbation theory for gauge equivalent equations. Sec. 5 
is devoted to the adiabatic approximation. Here we consider as an example the linear perturbation (excess 
gain or fiber loss) and the Raman self-frequency shift [19-21]. Concluding remarks are contained in the last 
section. 

2. Riemann-Hilbert problem for MNLSE 

In this section we summarize the basic results concerning the approach to the MNLSE based on the RH 
probem. Let us write the MNLSE in the general form 



iqz + ^q TT + ia {\q\ 2 q) T + P\q\ 2 q = 



(2) 
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where a and are real parameters. Eq. (2) is integrable via the 1ST method and can be represented as the 
compatibility condition U z — V T + [U, V] = for the following system of two linear matrix equations 

$ T = A(ft)[<7 3 , *] + 2ikQ<$> = 17$ - A(fc)$cr 3 , (3a) 

$ z = 0(fc)[a 3 , $] + [ — k 3 Q + 2ik 2 Q 2 a 3 - —HQ + kQ r a 3 - 2iakQ 3 ) $ 
\ a a ) 

= V$- fi(fc)$cr 3 . (36) 

Here the Hermitian matrix Q = ( _ ^ J represents the potential of the linear spectral problem (3a), 



q 

A(fc) = — (2i/a) (fc 2 - /?/4), fl(fc) = -(4i/a: 2 )(fc 2 -/3/4) 2 , the bar stands for complex conjugation and k is a 
spectral parameter. 

As seen from Eq. (3a), the MNLSE (2) belongs to the class of equations integrable by means of the 
Wadati-Konno-Ichikawa spectral problem [11]. However, Eq. (2) is not canonical among the equations of 
this class [16]. The canonical equation 

iq' z + \q' TT - iaq' 2 q' T + [3\q>\ 2 q> + a 2 |a'|V = (4) 

does not admit as obvious a physical interpretation as the MNLSE, but possesses the Lax representation 
too, 

$; = A(fc)[<7 3 , *'] + fakQ' + taQ' 2 ^ $' = U'& - A(fc)$'<r 3 , (5a) 

& z = n(k)[<j 3 , $'] + ( -k 3 Q' + 2ik 2 Q' 2 a 3 - % ^-kQ' + kQ' T a 3 + ^[Q' , Q' T ] 
\ a a 2 

+ l -^-Q' 4 ^ #' = v'& - si(k)&* 3 , (56) 

where Q' = ^ ^ ) The spectral problem (5a) associated with Eq. (4), as distinct from the spectral 

problem (3a), is compatible with the canonical normalization condition $'(fc = oo) = I, where I is the 
2x2 identity matrix. Eqs. (2) and (4) are gauge equivalent equations interrelated by the following gauge 
transformation: 

Q = g- 1 Q'g 1 (6) 

where g(r, z) — $'(fc = 0). Thereby, solutions of the MNLSE can be obtained from those of Eq. (4) by means 
of simple algebraic transformation (6). The RH problem formalism can be developed equivalently for either 
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the MNLSE (2) or Eq. (4), the RH problem data being invariant under the gauge transformation. Because 
the formulation of the RH problem with the canonical normalization condition has a number of technical 
advantages in calculation of soliton solutions, we will develop the RH formalism for Eq. (4) and give the 
transition relations to the MNLSE. 

To construct the RH problem associated with Eq. (4), consider the matrix Jost-typc solutions J'± of 
Eq. (5a) which satisfy the asymptotic conditions J± — > I at r — > ±oo. By the standard analysis of the 
Volterra-type integral equations for J' ± which follow from Eq. (5a) and the above asymptotic properties, we 
conclude that the following matrix function ((J±).; means the Z-th column of J' ± ) 

$;(*;) = (WiW, {j-uw), (7) 

being a solution of Eq. (5a), is holomorphic in the two quadrants of the complex fc-plane which are defined 
by the condition aRc(fc)Im(fc) < 0. The scattering matrix S'(k) is defined in the usual way: 

J'_E = J' + ES', E = exp (A(fc)a 3 r) . (8) 

Note that detS" = 1 due to det J'± = 1. The Zakharov-Shabat factorization [22] of the scattering matrix, 

allows us to represent <&' + {k) in two equivalent forms: 

$V = J'+ES'+E- 1 = J'_ES'_E-\ (10) 

Since the potential Q' is Hermitian, we have the following identities: 

J'i(k) = J'±\k), S'\k) = S'~ 1 (k), k G {k : Re(fc)Im(fc) = 0}. 

Hcrmiticity of the potential also enables to define the matrix function conjugated to and holomorphic 

in the rest two quadrants of the complex /c-plane, i.e., which are given by the condition aRe(k)lm(k) > 0: 

&_-\k) = & + \k) = ((j' + -\(k), (jl- 1 )^ , (ii) 

where (J'± 1 )/. denotes the l-th row of the matrix J' ± 1 and superscript t means transposition. The linear 
spectral problem (5a) possesses the parity symmetry [16]. It can be summarized to the following important 
identities: 

VJ' ± (k) = J' ± {k), VS'{k) = S'{k), VS' ± {k) = S' ± {k), 
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where V is the parity operator defined by VF(k) = F^(-k) - F {o ff\-k) = a 3 F{-k)a 3} F^ and F^fV 
are diagonal and off-diagonal parts of a matrix F. These identities give 

m*> = *<*>. *(*)-(£(*'; titXtiil) <12 » 

Here $'(fc) is a matrix function piecewise meromorphic in the complex fc-plane and discontinuous through 
the curve k e {k : Re(/c)Im(fc) = 0}. 

Let us return to the gauge equivalence between Eqs. (2) and (4). From Eq. (12) taken at k = it follows 
that the gauge transformation matrix g is diagonal, while Eqs. (5a) and (10) lead to the following expression 
for g: 



g = ^'(k = 0) = 



VI2 



V 



(13) 



exp^-ia / dr|g'| 2 ^ 

expl — ia J dr\q 

\ — oo 

The matrix function <!>(fc) which results from <&'(£;) by means of the transformation 

®(k)=g- 1 <S>'(k), (14) 

is a solution of the linear problem (3a), possesses the same conjugation and meromorphic properties as 
<£'(fc) does, and is an eigenfunction of the parity operator V too. For the linear problem (3a), the Jost-type 
solutions J±(k), the scattering matrix S(k), and its factorization S±(k) are constructed similarly to the 
problem (5a). The matrix function <I>(fc) is given through the Jost-type solutions and the factorization by 
the the same formulas as $'(fc) does (see Eq.(10)). The following relations are valid 



where g± = lim g. 

T^±00 



J±=9~ 1 J±9±, S = g+ 1 S'g-, S ± =g ± 1 S' ± , (15) 

Now we can formulate the RH problem associated with Eq. (4). Indeed, using Eqs. (7) and (11) as well 
as the relation E(k) = E{k) for k e {k : Re(fc)Im(fc) = 0}, we write (s' 12 (k) = S' 12 (k)) 

$'_ _1 (fc)$V(fc) = E(k)S' + \k)S' + (k)E~ 1 (k) 

= E(k) ( -^ {k) S[ f ] ) E~\k) ee G(k), (16a) 
$'(fc) -» I, k -» oo, (166) 
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where k £ {k : Re(fc)Im(fc) = 0}. It is a problem of analytic factorization of the nondegenerate matrix G{k) 
given on the contour defined by the following disconnected oriented set: 

C a = sgn(a) ({(oo, 0), (0, 0)} U {(0, 0), (-too, 0)} 



U{(-<x>,0),(0,0)}U{(0,0),(0,ioo)} 

of the fc-plane axes. Here sgn(a) means bypassing the contour in the reverse direction for a < 0. The 
functions 3>±(fc) are just a solution of the RH problem (16). The uniqueness of this solution is provided by 
the canonical normalization condition (16&). 

Substituting the asymptotic decomposition of <j>'(fc) in the inverse power series of k, 

$'(fc) =I + fc- 1 $' [1] + ..., 
into the spectral equation (5a) and taking advantage of Eq. (12), we reconstruct the potential Q': 

Q' = -[(73, = - lim k[a 3 , $'(fc)]. (17) 

a a k— »oo 

Hence, to solve the MNLSE (2), we should at first solve the RH problem (16), then obtain the potential 
Q' (17) and transform it by Eq. (6). It should be noted that the RH problem (16) remains unchanged 
under the gauge transformation (14), except for the normalization of the matrix $(fc) at infinity, namely, 
—> 9~ l at k — > oo. 

In general, the functions det$' + (k) and det$>'_~ (k) have zeros in their regions of analiticity, the RH 
problem being said to be nonregular (or with zeros). It follows from Eq. (11) that zeros of the above 
determinants are complex conjugate, while the parity symmetry (Eq. (12)) tells us that zeros appear by 
pairs, i.e., det$' + (±fcj) = for the j-th zero kj. The solution of the RH problem with zeros can be 
factorized [23] 

*' ± (k) = *' 0± (k)T(k), det$' o± (fc)^0 (18) 
by means of the solution $' ± of the regular RH problem: 

(K-r\k)K + ( k ) - r(fc)G(fc)r- 1 (fc), (i9«) 

$' (fc) -» I, k -» oo, (196) 

which is posed on the same contour C a . The matrix function T(k) represents the contribution of zeros. In 
the case of the single pair of zeros ±k\ = fc±i this function is given by (see also Ref. 18) 

j,l=±l k ~ kl j,l=±l 3 
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where D ln = (k n - k t ) 1 (p l \ p n ), D.. 1 = (D x ).., (p t | p n ) = (p;)i(Pn)i + (Pi) 2 (Pnh, vector-rows (p ±1 | are 
related to vector-columns | p±i) by conjugation, i.e., (p±i \ — \ p±i)^ , and the latter are given by 

*'+(*±i) \P±i) = 0. (21) 

The parity symmetry immediately gives a relation between | p_i) and | Indeed, the identity & + {ki) \ 

p+i) = cr 3 <i>'(A:_ 1 )(T3 | p+i) = leads, due to uniqueness of | P-i), to the relation | = 03 \ p+i). From 
Eq. (11) it follows that (p±i \ $>'_~ (fc±i) = 0. This identity and Eq. (21) ensure that the matrix functions 
$' 0+ (fc) and ($Q_) _1 (fc) are holomorphic in the respective quadrants of the complex fc-plane. 

In the general case of the RH problem with N pairs of zeros of det$>' + (k) the matrix function T(k) is 
given by Eq. (20) with j, I e {-N, -N + 1, -1, 1, ...N - 1, TV} and fc_ s = -fc s . 

3. Soliton solution of MNLSE 

The RH problem data are divided into two parts: discrete data {kj, 
| p+j),j = 1, --,N} (2N is the whole number of zeros of det$' + (k)) and continuous datum G(k). Soliton 
solutions correspond to the RH problem with zeros provided G(k) = I, i.e., & ±{k) = I. In other words, only 
the matrix T(fc) is responsible for solitons. We will consider the simplest case of one-soliton solution of the 
MNLSE (2). We have from Eqs. (6) and (17): 

q = -r- x (fc = o) lim % 3 , r(fc)]r(fc = o), (22) 

a fe^oo 

where the matrices r(fc) and r _1 (A:) are given by Eq. (18) and it is taken into account that g = T(k = 0). As 
regards the coordinate dependence of the vector-column | p+i), it is determined by differentiation of Eq. (21) 
with respect to r and z and taking advantage of Eqs. (5). It should be noted that Eq. (21) determines | p+i) 
only up to an arbitrary norm. We obtain in a particular case 

| p +l ) T = A(fci)<7 3 | p+i), | p+i) z = fi(fci)a 3 I P+i>- (23) 
Integration of the above equations gives 

| p+i) = e /CT3 | p°), \P°) = (fo), 

where / = A(fci)r + Q(ki)z and | p°) is an integration constant determined up to an arbitrary norm. Let us 
define Pi/p% = exp(a + iip), where a and ip are real constants. Ultimately, we have 

ip«>-( cxp( " e + -r +/) ). /=-!(*?-i>-S(*?-i) 2 - ^ 
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Substituting Eq. (24) into Eq. (20), we obtain for T(fc) the following expression 



r(fc) = i- 



2e° 



he* {D~\ 



he** (Dll + Dz\) 



(25) 



Here 



8£r? / 4 
x = a + 2Ref = a - — r + - 
a \ a 



e-v 2 - 



(3 



(26) 



ip = if + 21m/ = if — 




tl 2 



4£ V I z 



and ki = £ — it] (due to the condition aRe(fci)Im(fci) < we will have £ > and r\ > for a > 0). From 
the definition Di n = (k n — ki)^ 1 {pi \ p n ) we derive the following properties of the matrix D: D ++ = — D—, 
D_ + = — D + -. Hence, 

D+l ±DZ\ = - (D ++ ± D_+) (dcti?)- 1 = (D ++ T D^+y 1 



= e a [ — chx =F -rsha; 



(27) 



Now, we have all to calculate one-soliton solution q s of the MNLSE (2). Substituting Eqs. (26) and (27) into 
Eq. (25), then, in its turn, Eq. (25) into Eq. (22) we obtain 



Qi 



a (kie x + fcie-*) 2 



(28) 



This is a general form of the one-soliton solutions to the MNLSE (2) which depends on four real parameters 
£, r], a, and ip. The solution (28) is written in terms of the coordinates x and ip (26) comprising linear 
combinations of the normalized retarded time r and distance along the fiber z. It should be stressed that 
a enters the denominator of the soliton solution (28) . In other words, we account nonperturbatively for the 
pusle self-steepening effect. Similarly to the NLSE, the parameters a and ip play the role of the initial position 
and phase, respectively, while the other two parameters £ and r\ do not admit so obvious interpretation. In 
any case, we see from Eq. (28) that the normalized half-width w and velocity v of the soliton are represented 
by 

w = -£-, v=P--(e-V 2 )- (29a) 



a 



a 
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As regards the soliton amplitude A, it is natural to admit it in the following form: 

A = — — -. (296) 

The MNLS soliton (28) has a number of peculiarities which distinct it form the NLS soliton. First, the 
MNLS soliton has nonzero phase difference at its limits. Indeed, 

k 1 e~ x + kie x I ki/kf, x—>oo, 



7-1 



[k x e x + kxe- x ) 2 [ k 1 /k 1 , x -oo 

and the said phase difference reads 

Alp = arg (q s (z —> oo)) - arg (q s (z -oo)) = 6arg(fci). 
Further, the integral of the soliton amplitude 

n=0 



oc 

/ 



where (j) n = (j + l) ( j + 2) • ... • ( j + n — l), depends on the parameters £ and r/, in contrast to the same 
integral of the NLS soliton which does not depend on the soliton parameters. Moreover, the important 
invariant of Eq. (2), namely, the number of particles or the optical energy of the soliton, 

oo 

E= [ dr\q s \ 2 = -axg(Ai), < |arg(*i)| < J, 
J a 2 

— oo 

has the upper limit 27r/|a|. The phase difference and optical energy of the MNLS soliton are related: 
= -3Ea/2. 

The above properties of the MNLS soliton (28) resemble those of the dark NLS soliton, which also has 
nonzero phase difference and a relation between its energy and the phase difference [24]. Nevertheless, the 
soliton (28) with (3 > reduces to a bright NLS soliton at a — > 0. To carry out this limit one should take into 
account that the Lax pair for the NLSE should be produced at a — ► from the Lax pair (3) for the MNLSE. 
This condition implies that the spectral parameter k depends on a and gives the following prescription: 
2 ( k liNLS " P/^)/ a -k NLS at a -» or 



(< 2 -t) 



-t ^ - » 

SOI '/07 



a a 
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where Izmnls = £ — if] an d fcjvLS = Co + *?7o ln other words, we have the decomposition 

e= ^2^° + ° (a2) ' ^2^ ?7o + 0(a2) ' (30) 

which transforms the MNLS soliton (28) to the NLS soliton 

q s = (2i?7 /V/3)e l,/ ' o sechx , where ip Q = ip + 2£ Q r - 2(g - Xq = a _ 2 ??0 (t - 2£ z). 

4. Perturbation-induced evolution of RH problem data 

Following our strategy, to find correction to the soliton solution of the MNLSE caused by a perturbation, 
we should at first derive the perturbation-induced evolution equations for the RH problem data. Decompo- 
sition of these equations in the asymptotic power series with respect to the perturbation will produce the 
consequent corrections to the soliton solution. 

In Sec. II we have seen that the gauge equivalent Eqs. (2) and (4) have resembling 1ST schemes with 
simple mutual relations (15). As the RH problem data are invariant under the gauge transformation (6), 
we can choose between the two 1ST formulations the most convenient one for calculation of corrections to 
soliton solutions. Consider eq. (2) with a small perturbation at the r.h.s. of the equation, i.e., Eq. (1). 
The perturbation causes a variation 5U of the potential, what, in its turn, leads to a variation 5 J± of the 
Jost-type solutions. From the spectral problem (3a) we obtain the following equation for 6 J± 

(SJ±) r = -A(k)6J±a 3 + 5UJ± + USJ±, lim (J^S^) = 0. 

Therefore, 



5J± = J±E y j &tE- 1 J ± 1 5UJ±E^ E-\ 



(31) 



where SU — (SU/Sz) 5z. It should be stressed that the same representation of the variation SJ' ± of the 
Jost-type solutions to the spectral problem (5a) follows from Eq. (31) by the simple substitutions U — > U' 
and J± — ► J'± and such a procedure can be carried out at any step below. 
Now wc introduce a useful matrix function 



T 

/err 
AtE- 1 ^ 1 — 



(32) 



Then the variation derivative 5J±/5z takes the form 



^± = J±ES±j(±oo,T)Sg 1 E- 1 . (33) 
oz 
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From the evident relation S = lim E 1 J-E we obtain the variation of the scattering matrix: 



^=S +7 (-oo,oo)Sl\ (34) 



as well as of its factorization (9): 



" H =^+7(-oo,oo)M 22 , M 22 =(l J). (35) 



Sz + ,v ^ ^ o 1 

What regards the variation 5<fr + (k) for k G {Rc(fc)Im(fc) = 0}, we obtain from Eqs. (10), (33), and (35): 

= ^ ( J +^ 5 + S_1 ) = T ) + °°) M 2 2 ) ^ 

= $ + £^-7(t,oo)Mh +7(-oo,t)M 22 ^" 1 , 

where Mn = diag(l,0). Hence, in the case of perturbation the evolution equation for $ + with respect 
to z (3a) gains an additional term responsible for the perturbation: 

{^ + ) z = V^+-n{k)^ + <T 3 + ^ + EUE~ 1 , (36) 

where the fundamental matrix function II(fc) is defined as 

n W - +7( -„,,)m 22 . ( ) . (37) 

Eq. (36) was derived for k belonging to the curve {Re(fc)Im(/c)} = 0, involved in the formulation of the 
RH problem, but it can be analytically continued into the domain {k : aRe(fc)Im(fc) < 0}, where is 
holomorphic. Indeed, from the linear problem (3a) we get 

($; 1 (5$ + ) t = $- 1 5U<f>+ + A(k)[a 3 , 

Hence, the integrand of Eq. (32) is represented as (E~ 1 $+ 1 (5$ + /5z)E) t and the matrix function 7(±oo,r) 
reads 

7(±oo,r) = E-^^±E ~S^, (38) 
oz oz 

where we have used Eq. (10) to compute the limits r — ► ±oo. Now, in virtue of the explicit structure of the 
matrices S± (9) and analiticity of their diagonal elements in the domain {k : aRe(k)lm(k) < 0}, we conclude 
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that the entries of the matrix H(k) (37) are functions meromorphic in this domain and having simple poles 
at the zeros of the RH problem, i.e., at the zeros of dct$ + (fc). Analiticity of $+(fc) in the above domain 
gives immmcdiatcly an important identity 

^+(±k j )E(±k j )Res{U(k), ±kj} = 0, (39) 

where Res{-, ±fcj} stands for the residue at k = ±fcj, det$ + (±fc :) ) = 0. 

In order to derive the perturbation-induced evolution of the RH problem data, let us consider for simplicity 
the case of a single pair of zeros ±k\. The full set of the independent data involves in this case k\, \ p+i), 
and G(k). By the definition of ki, we have 



o = A det $ +(fcl) = | dct<Mfcl) + ^ ^ det$+ 



k—ki 



trn(fc)dct$+(fc)) + ( ^-det$_ 

Jk=k 1 dz \ dk /k=k 1 

Here we took advantage of Eq. (36) and of the identity trV = 0. Decomposition of n(fc) into the regular and 
singular parts at k = k\, 

n( fc) = n r ( fc) + Res{ fc n( ^ fcl} , (40) 

and use of the formula det$ + = dct$ 0+ dctr = dct<I> 0+ (fc 2 — k 2 )(k 2 — k^ 1 which follows form Eqs. (14), 
(18), and (20) with <f> = </ _1< 3?o help us to compute the limit k — > k\: 

dfci _ / tr [n r (fc) + (fc-fc 1 )- 1 Res{n(fc),fci}] (P-kDiP-kl^dct^o+iki 



dz \ d_ 

8k 



(k 2 -k 2 )(k 2 -/ci)- 1 det$ 0+ (fc) 



k—k\ 

= -Res{triT(fc),fci}. (41) 

The perturbation-induced evolution equation for the vector | p +1 ) can be obtained by differentiation of 
Eq. (21): 

d -* + (k)) |p +1 )+$ + (fc 1 )A| p+l)=0 . (42) 



dz J k=k 1 dz 

Here the full derivative of $+(fc) with respect to z is given, due to Eq. (36), as follows: 
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-n(k)<f>+{k)a 3 +<f>+{k)E(k) [IL r (k) 



dfci / d , ,, . 



Res{n(fc),fci} 
k — k\ 

I P+i>, 



^(fc) 



P+i) 



k—ki 



(43) 



k—k\ 



where we have used the identity (21) to exclude an evidently vanishing term: V{k 1 )<^ + {k 1 ) \ p +1 ) = 0. 
Turning again to Eqs. (36) and (40), we obtain the following relation: 

[(k - fci^^XS+Mfc)]^! P+i) = E(k 1 )Res{U(k),k 1 }E- 1 (k 1 ) | p +1 >. 
On the other hand, the l.h.s. of this relation can be transformed in that way: 



[{k-k 1 )^ 1 (k){9+Uk)] k=ki \p +1 ) 



dfci 
dz 



P+i) 



dJfei 



- ([(fc - fci) S+^fc)] J fc=fci $+(fei) | p+i) = | P+1 ), 
where we have once again used Eq. (21). Hence, we have derived an important identity 

£7(fci)Res{II(fc), fci}^- 1 ^!) | p+i) = — ^ | p+i). 



Using now Eq. (39) we can write 



Z + (k)E(k) Res{U{k }> kl} E-Hk) 
k — k\ 



P+i) 



k—ki 



d_ 

dk 



J B(fc 1 )Rcs{n(fc), fcr}^- 1 ^!) | p+i). 



k—k\ 



Collecting Eqs. (42)- (45), we obtain the perturbation- induced evolution of the vector | p+i), 

A | p+1 ) = o( fcl )a 3 1 p+i) - s^on^fci)^- 1 ^) I p +1 ), 

or, for r-independent vector | p°) = F _1 (fci).E _1 (fci) | p+i), where F(k) = exp I J dzVl{k)a 3 



(44) 



(45) 



— \p°) =-F-\k 1 )Il r (k 1 )F(k 1 ) |p°). 



(46) 
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Finally, we should derive an evolution equation for the continuous datum G{k). Eqs. (36), (11), (16) (written 
for the quantities without prime), and (36) give 

G z {k) = n{k)[^,G{k)]+G{k)E{k)I{{k)E- 1 {k)+E{k)TlHk)E- 1 {k)G{k). 

This equation is simplified by introducing the r-indcpendent matrix G° = F~ 1 E~ 1 GEF: 

G°{k) = G°(fc)F- 1 (fc)n(fc)F(fc) + F- 1 (fc)n t (I)F(fc)G°(fc). (47) 

It is evident that T-independence of the l.h.s. of Eqs. (41), (46), and (47) enables us to simplify further these 
equations by taking one of the limits r — > ±oo, at which the fundamental matrix II(fc) (37) has only one 
nonzero column. 

As was pointed above, we can use the gauge equivalent 1ST schemes for the derivation of perturbation- 
induced evolution equations for the RH problem data, due to the gauge invariance of the latter. Indeed, 
from Eqs. (38), (14), and (15), using the parity symmetry, one can easily find the folowing relation between 
the gauge equivalent fundamental matrices: 

n'(fc) = n(fc) - ^- 1 (fc)$; 1 (fc)n^(fc -» oo)$ +{k)E{k) 

= U(k) + E-\k) K \k) ( ^°°) _ 722 (°oo,.) ) k ^J >+(k)m - 

where II' (k) is associated with the canonical equation (4). It is easy to see that the second term in the r.h.s. 
of this expression for H'(k) gives no contribution to the evolution equations (41), (46), and (47). For instance, 
Rcs{trn'(fc), ki} — Res{trII(fc), k\}. This once again ensures that one can use either of the gauge equivalent 
fundamental matrices II or II' expressed through <f> + , SU/Sz or $' + , 6U'/5z, respectively, for calculation of 
corrections. It will be more convenient to use the former fundamental matrix II, because precisely 5U\2/8z 
is proportional to the r.h.s. of Eq. (1), treated as the perturbation to the MNLSE. 

5. Adiabatic approximation 

It should be stressed that the evolution equations (41), (46), and (47) of the preceding section arc 
exact, but just this circumstance prevents their direct use for concrete calculations. They include unknown 
function <& + (/c), a solution of the perturbed equation (36). The way to determine $+(fc) explicitly defines 
the sort of the approximation used. Here we develop the simplest variant to account for perturbation, the so 
called adiabatic approximation. In the adiabatic approximation, soliton shape is considered to be unchanged 
instantaneously under the action of a small perturbation, while the soliton parameters, being constant in 
the integrable case, acquire slow z-dependence. In other words, we impose the condition G(k) = I and, 
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consequently, & +( k ) = 1 ( 19a )- Thcn > Ec l s - ( 14 )> ( 13 ) and ( 18 ) S ivc $ +( fc ) = .9~ 1$ +( fc ) = r_1 ( fc = O)r(fc). 
Let us define the perturbation matrix R as 

^ = 2ik^ = 2kR. (48) 
oz oz 

Then the r.h.s. of Eq. (1) is nothing but Ri 2 . As was mentioned above, taking the limit r — > ±oo considerably 
simplifies the evolution equations (41) and (46) for soliton parameters. From Eq. (32) we obtain 

oo 

7(fc) = 7(-oo,oo) = 2k J dTE- 1 (k)T-l(k)T(k = 0)RT-l(k = 0)T(k)E(k). (49) 

— CO 

At the limit r — * oo the fundamental matrix U(k) defined by Eq. (37) takes the form 

»-•<*>- Mil L.- <») 

With this representation for II(fc), the evolution equation (41) is written as dfci/dr = — Rcs{722(fc), fci}. 
Taking into account Eqs. (25) and (27), from Eqs. (41) and (49) we obtain 



dz J (kie~ x + k^) 



iakf / dxe x " w _ v (51) 



where r (x, z) = e "M^V^ z) and the variables a; and ip are given by Eqs. (26) with account of possible 
z-dependence of the soliton parameters £ and r/: 



8HV 32 / / 2 2 /A 



2 



2„2 



4^ 



For the soliton parameters £ and 77, fci = £ — ir/, we get from Eq. (51): 



00 , — 2 
d£ za /" fc^e 21 — fc x e 2 



ia f , fcf 
= — / da; 



dz 2 j/ (fcie- x + fcie x ) 



(r (x)+r (-x)), (52) 



15 



dr/ 



a f , k 2 e x 
= — — / dx 

2 i (fcie" : 



—2 



e 21 + A^e 



(fcie _a: + kie x )' 



(r a (x) +r Q (-x)). 



(53) 



It is noteworthy that only the symmetric part of the perturbation, i.e., such that r Q (— x,z) = r a (x,z), 
changes the soliton parameters £ and rj. 

In order to obtain the adiabatic evolution for the rest two soliton parameters a and ip, we turn to Eq. (46). 
Rewriting it by entries with the matrix II given by Eq. (50) and taking into account Eq. (40), we obtain 



dz 



z 

(X|> ( 2 / : )~<i-.:/'l ;/'!•• -^r = - : : l>2 



and for p\lp° 2 = e a +^: 



dz 



c — e Jr22 



z 

(fci)-exp^-2 J dzfi(fc 1 )^7 rl2 (fc 1 ). 



(54) 



Now we should extract the regular part of 7(fc) according to the decomposition 7.2 (k) = 7 r . 2 (/c)+Res{7.2(fc), fei}(fc- 
fci) -1 (40). After long but not difficult calculations with account of the explicit expression for T{k) (25) we 

get [C = a- (32/a 2 ) / dz ^ (£ 2 - V 2 - 0/4) 



z 00 

7n 2 (fci) = ^expf2 dz^k^e^ f ^ 

\ J J J (fcie-* + fci 



(2k 2 C - 3%rj) (r (x) + r (-x)) + (£ 2 - rf - 2k 2 x) (r (x) - r (-x)) 
+2|fc 1 | 2 r (x)e- a; |, 



7 r22 (fci) = 7^ / < l < 



r p (x) + r Q (-x) 



(|fci| 2 e- x -3^r?e x ). 



Now from Eq. (54) we obtain 



a 



da 

dz" 2€n J ( k 

— no V 



00 

/ K Zk x J C (^ + tie ~ X ) {r ° 

J (kie x + kie x ) I x 7 



(x) + r Q (-x)) 



16 



dip 



i£v (he x + he x ) + x [k\e x + k x e~ x 
dx 



(r (x) - r c (-a;)) j, 
C (k\e x - k\e~ x ^ (r a (x) +r (-x)) 



(55) 



dz 2 ^ J (k ie - x + k ie x ) 

— oo ' 

kie x (k\ - + k x e~ x (k\ + i£r]) - x (k\e x - ~k\e~ x ^j 



(r a (x) - r (-x)) 



(56) 



For the case of the symmetric perturbation, i.e., when r (— x, z) = r a (x,z) holds, the evolution of a and ip 
is determined by the evolution of k\: 



da 



C 



dz 



dz 



dp C . 



— = -— Im<^ ki— - \ , -p- = — Reiki— ^- > ■ 



dz ^77 



dz 



(57) 



Hence, in the symmetric case the adiabatic evolution of the soliton parameters is determined by a single 
integral involved in the r.h.s. of Eq. (51). It is remarkable that the symmetry condition r (— x, z) — r Q (x, z) 
is rather natural for the MNLSE, because the soliton (28) itself satisfies this condition: q s (—x, z)e~ 1 ^ — 
q s (x,z)e-^. 



Example 

As an example of using the above formalism, we analyse here in the adiabatic approximation the effect of 
two different perturbations to the soliton propagation. We will consider a perturbation comprising two terms, 
namely, the linear perturbation which describes the fiber loss or excess linear gain and the Raman driving 
which is reponsible for the soliton self-frequency shift effect [19-21]. Such a choice of perturbation types 
does not in any way exhaust the set of all perturbations relevant for the femtosecond soliton propagation in 
optical fibers and serves mainly as an illustration of the perturbation theory developed above . The results 
will be presented in a form suitable for numerical calculations. 

Let us write the perturbation in the following form: 

r = ieiq + s 2 \q\lq, (58) 

where e\ and £2 are small parameters characterizing the linear gain (ei > 0) or the fiber losses (ei < 0) 
and the Raman driving, respectively. Substituting the soliton solution (28) into Eq. (58), we obtain the 
characteristic perturbation r Q in the form 

fr> k x e~ x + kie x 

r Q = -4ei — = — 

a (k-ie x + k x e- x Y 
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+ie 2 



-2x 



(59) 



\aj (kie x + k ie - x Y (k x e- x + k x e x ) ' 

Note that this perturbation is symmetric in the sense of the preceding section, i.e., r Q (— x, z) = r {x,z). 
Hence, evolution of the soliton parameters a and ip is determined by evolution of k\ in accordance with 
Eq. (57). With the perturbation (59) Eq. (51) takes surprizingly simple form: 



dfci 
~dz~ 



k, 

-Ely (1 - 20cot(20) - 2i0) 



-16e 2 



fcilfcil 4 



1 - 20cot(20) - -sin 2 (20) 



ith 



d0 



= 2ei0. 



Here = arctan(r?/£), < < tt/2, cot(20) 



dz 

-V 2 )/(^V), and sin(20) 



(60) 



(61) 



It should be noted that the parameter a enters Eq. (60) in a denominator representing thereby a non- 
perturbative account of selfsteepening of the soliton. In the limit (30) Eq. (60) reduces to the well-known 
equations [5] for the NLS soliton parameters. From Eq. (60) one can derive evolution equations for the 
soliton parameters (29): 



dz 



%(l-20cot(20)) + 32 £2 ^^ 
2 or 



l-20cot(20)-isin 2 (20) 



dv 
d~z 



ei(l-20cot(20)) +32£ 2 



fell fci| 4 



l-20cot(20)--sin 2 (20) 



V 

a 



-8ei0u;" 



dA 

d7 



(60cot(20)-l) - 16e 2 



fci I 4 / 1 

' 1 -20cot(20)--sin 2 (20) 



a 



A. 



It follows from Eqs. (61) that the soliton energy E = J dr\q\ 2 = 40/|a| does not depend on the Raman 

— oo 

driving: 

— = 2 £l E, 
dz 
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just as in the case of the NLSE [5]. 
6. Conclusion 

We have developed a perturbation theory for the MNLS soliton. It has been shown recently [25] that the 
MNLSE utilizing the notion of slowly varying envelope is still valid up to 3 - 5 periods of field oscillations 
within the envelope. Our formalism is based on the RH problem associated with the Wadati-Konno-Ichikawa 
linear spectral problem. We have demonstrated previously the efficiency of the similar approach to the 
perturbed Manakov system [15] , as well as to nonlinear evolution equations integrable by the NxN Zakharov- 
Shabat spectral problem [14]. It is not mere a chance that the formalisms for both linear and quadratic 
spectral bundles resemble each other. Indeed, the main idea behind the RH-based approach consists in the 
analysis of analytical behavior of the Jost-type solutions to the spectral problem. Besides, we have stressed 
the gauge invariance of the RH problem data. Just this property is responsible for the definite unification 
of the objects used in both cases. In the recent paper [26] the RH problem has been used for studying 
asymptotics of the MNLSE solution associated with the continuous spectral data. 

As regards the examples of perturbations considered above, we confine ourselves to the adiabatic approx- 
imation. This restriction is by no means concerned with the present approach. Moreover, we have derived 
the perturbation-induced evolution equation (47) for the continuous RH problem datum. It is precisely this 
equation which is necessary for the description of the soliton shape distortion and emission of linear waves. 
This problem will be solved in a forthcoming paper. 
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